Stop doing full-array work at open to answer approximate questions - #23

Stacked on #22 — review that first, and this diff is just the two commits on top.

Opening an NWBSlicer did several full passes over a recording's timestamps to compute statistics that never needed that precision. On a 5-stream 5.7 GB recording:

The cold number is the one that matters for batch analysis, where every file is seen for the first time.

What was slow

np.median over every inter-sample interval — 45 ms a go, several times per stream. infer_nominal_rate uses the median only to centre a trim window; the estimate it returns is a mean over the intervals that survive. Partitioning 7 million floats to place a window is precision nobody asked for. clockmodel's _nominal_period and _auto_gap_threshold do the same thing on the same arrays.

np.var computed twice to evaluate one or. Short-circuiting meant the second full pass only ran when the first comparison failed — which is exactly the irregular-stream case, i.e. the slow one.

electrodes.table.to_dataframe() — 12 ms per electrical series. It materialises every column the writer stored (position, group, filtering, …) into pandas, and resolves an object reference per row, in order to read one column.

Why subsampling the median is exact here

Timestamps arrive on a device's quantisation grid, so their differences cluster onto that grid and the median snaps to a grid point that a stride cannot move. On a real 7.09M-interval stream, strides of 16 through 1024 all return 3.3301000002e-05, and infer_nominal_rate returns 30147.9172 from every one of them.

Verified end to end on that recording: all five streams reconstruct bit-identical dejittered timestamps, and the gap threshold matches to the last digit. Two of those streams reach the inference path with non-trivial rates (30137.773554857125 and 30000.082189805245), so this isn't only exercising the easy cases.

The helper is named interval_median and lives in util.py, which imports nothing local — clockmodel can't import from slicer, since slicer imports clockmodel. The name states what licenses the shortcut: intervals inherit a grid, and other quantities don't.

What deliberately did not change

_auto_gap_threshold keeps its percentile on the full array. Subsampling the whole function is 16x faster and moves the threshold 1.4% (0.013372 → 0.013557), because a 99.9th percentile over 200k values rests on their top 200 instead of on 7,000. That threshold decides which interval jumps count as real gaps, so moving it re-segments the stream and re-times its samples — a behaviour change wearing a speedup's clothes. A test asserts the asymmetry so nobody finishes the job later.

through_c's median is untouched for the same reason: it averages fit residuals, which are continuous and have no grid to snap to.

Tests

tests/test_slicer_perf.py pins equivalence, not speed — interval_median against np.median across the threshold, _nominal_period and _auto_gap_threshold against themselves computed with a full median, and inferred rates unchanged. Because the trim is the only thing that makes a mean usable here, one test gives a stream 400 dropped-packet holes and asserts it still reports 30 kHz rather than an average including the holes. Another pins that the label column read honours the electrodes region, since a series may reference a subset of the table.

The fallback for a table with no label column now names channels from the id column rather than a DataFrame index. That is not a change — the index to_dataframe() built was the id column — but it has to be explicit now that the DataFrame is gone, because positional indices coincide with ids only when a series references the whole table in order.

256 tests pass; ruff clean.

What's left, for the record

After this, opening is ~119 ms, of which 108 ms is inside pynwb.read() building Python objects. About 24 ms of that is pynwb#2253 (submitted upstream): FeatureExtraction.__init__ copies times into a Python list, reading the backing dataset one element at a time — 23,630 separate h5py reads on this file, 84x slower than one bulk read. The rest is hdmf's object construction and not ours to fix cheaply.

Chunk iteration was profiled separately and needs no work: per-chunk CPU in our code is ~0.07 ms against an 11–16 ms chunk, and raw h5py forward reads run at the same 1.09 GB/s the iterator achieves. It is disk-bound, not us.

🤖 Generated with Claude Code